Quantum phase transition in Bose-Holstein model in two dimensions 
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We derive an effective d-dimensional Hamiltonian for a system of hard-core-bosons coupled to op- 
tical phonons in a lattice. Away from half-filling, we show that the presence of next-nearest-neighbor 
hopping in the effective Hamiltonian leads to a superfiuid-to-supersolid transition at intermediate 
boson-phonon (b-p) couplings, while at strong-couplings the system phase separates. However, at 
half- filling and at a critical b-p coupling (as in the xxz- model), the system undergoes a superfluid- 
to-charge-density-wave transition without any signature of supersolidity. Our analyses is based on 
extensive calculations of the structure factor, the superfluid fraction, the Bose-Einstein condensate 
fraction, and the system energy at various fillings. We present a phase diagram for this system and 
compare it to that of the xxz-model. We also demonstrate explicitly that the next-nearest-neighbor 
hopping (in the absence of nearest-neighbor hopping) in the effective Hamiltonian leads only to a 
single transition - a first-order superfluid-to-supersolid transition. 
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I. INTRODUCTION 



The successful mimicking of an actual lattice using op- 
tical standing waves marks one of the most significant 
scientific advances of this decade^. The biggest advan- 
tage of this kind of an optical lattice is that the ratio 
of the kinetic energy and the interaction energy of the 
particles can be controlled at will. This has led to a 
flurry of activities among atomic and condensed matter 
physicists across the world. The excitement among the 
condensed matter physicists stems from the fact that it 
not only gives a testing ground for some of the most 
intriguing phenomena of nature predicted earlier, but 
it also paves the way for the discovery of new physical 
phenomena. For example, soon after the creation of a 
two dimensional (2D) optical lattice, it has been used to 
experimentally verify^ the predicted transition 3 - from a 
superfluid state to a Mott insulating state of a bosonic 
system. Another possibility is the verification of the the- 
oretically predicted supersolidity due to vacancies^^. A 
signature of supersolidity is the simultaneous presence of 
both diagonal long range order (DLRO) and off diago- 
nal long range order (ODLRO) 7-8 . There have not been 
many studies of this interesting phase of matter until it 
was recently observed in helium- l 9 i 10 . This discovery led 
theorists to study bosonic models in different kinds of 
lattice structure d 11 ' 12 ! 13 , and with various types of inter- 
actions among these particles^. 

There has been very little attention given to a system 
of bosons interacting with phonons. Recently Pupillo et 
al*£- have studied such a possibility where the bosons 
could be coupled to the acoustic phonons generated by 
polar molecules trapped to form a lattice. In this pa- 
per, we have considered a Bose-Holstein model compris- 
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ing of hard-core-bosons (hcb) coupled to optical phonons 
generated by the vibrations of the underlying lattice. 
An example of such hcb is a collection of tightly-bound 
Cooper pairs originating from electronic polarization 
processe o 16 i 17 . Additionally, strong-coupling between 
electrons and intcrmolecular-phonons also produces hcb; 
when such hcb couple to intra- molecular phonons, the 
system can be studied by a Bose-Holstein models. Start- 
ing with a minimalistic model, involving momentum in- 
dependent b-p coupling, we have derived an effective d- 
dimensional Hamiltonian for hcb by using a transparent 
non-perturbative technique. The region of validity of our 
effective Hamiltonian is governed by the small parame- 
ter ratio of the adiabaticity t/u>o and the b-p coupling g. 
The most interesting feature of this effective Hamiltonian 
is that, besides a nearest- neighbor (NN) hopping, it con- 
sists of next-nearest-neighbor (NNN) hopping and NN 
repulsion. Our approach gives a microscopic justification 
for the origin of these important additional terms. We 
study our derived effective Hamiltonian in 2D by using 
exact diagonalization technique. For exact diagonaliza- 
tion, we have used a modified Lanczos algorithm^ on 
lattice clusters with 4x4, \/l8 x y/18, and x V%0 
sites. We have shown that, except for the extreme anti- 
adiabatic limit, the hcb coupled with optical phonons 
can show supersolidity above a critical value of the b-p 
coupling strength. 

The paper has been arranged as follows. In Sec. HH wc 
have derived the effective Hamiltonian for a system of hcb 
coupled to optical phonons. We discuss briefly the ba- 
sic difference of this effective Hamiltonian with that for 
fermions 2 * 3 . Next, we apply mean field analysis to this 
Hamiltonian in Sec. IIIII and obtain a mean field phase 
diagram. In Sec. IIV1 we discuss in detail the DLRO by 
studying the structure factor. Here, we also present key 
numerical results. In Sec. [Vj we discuss two important 
quantities - the Bose condensate fraction and the super- 
fluid fraction. Sec. IVII deals with calculating the free 
energy of the system for different situations and parame- 
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ter values. The curvature of the free-energy- versus-filling 
curves is used in deciding whether the system phase sep- 
arates or not. Finally, in Sec. IVIIi we present the results. 



II. EFFECTIVE HAMILTONIAN 

We start with a system of spinless hcb coupled with 
optical phonons on a square lattice. This system is de- 
scribed by a Bose-Holstein HamiltonianSI 

H hoi = -ty^ y b] +s bj + up qtfflj + gu) o y^^nj(aj + oj), 

i,<5 3 3 

(1) 

where 5 corresponds to nearest-neighbors, Uj = bjbj with 
bj being the destruction operator for hcb (and not of 
electrons as in the Holstein model) , while (as in the Hol- 
stein case) Oj is the destruction operator for phonons, 
and ujq is the single vibrational frequency for simple har- 
monic oscillators. Then we perform the Lang-Firsov 
(LF) transformation 22 on this Hamiltonian which pro- 
duces displaced simple harmonic oscillators and dresses 
the hopping particles with phonons. It is important to 
note that although we are dealing with particles different 
from fermions, we can still perform the same LF trans- 
formation. This is because, under the LF transformation 

given by e s i?hoie _,s with S = —gJ2i n i( a i ~ a \)i bj and 
a,j transform (like fermions and phonons in the Holstein 
model) as follows: 

bj = e s b je - s = b j e- g ^- a i ) , 

hj = e s a,je~ s = a,j — grij. (2) 

This is due to the unique (anti-) commutation properties 
of hcb given by 

[Mi] = [h,b]]=0, iovijtj, 
{b h b\} = 1. (3) 

Next, we take the unperturbed Hamiltonian to be given 
by2i 

H = up ^ QjQj - 9 2uJ o^2 b ] b j - J i^2( b 3 b 3+s +H.c), 

3 3 3 

(4) 

and the perturbation to be 
H' = J2 Hj = -Ji Y^( b h + s{s(si 1} + H.c), (5) 

3 3 

where S± = exp[±g(aj — a,j + s)],Ji — texp(— g 2 ), and 
g 2 oJo is the polaronic binding energy. Here, Hq + H' con- 
stitutes the LF transformed Bose-Holstein Hamiltonian. 
We then follow the same steps as in Ref. [20] to get the 
following effective Hamiltonian in d-dimensions for our 




FIG. 1: Comparison of J%jJ\ and J z /Ji at two different val- 
ues of t/uo and various values of g. 

Bose-Holstein model 

H eff = -9 2 u n 3 - J i b ) b 3+ s 

3 3,S 

-J 2 Y b )+S' b 3+s - 0.5J z ^2rij(l - rij +s ), 

j,S,S'^6 j,5 

(6) 

where J z = ( J 2 / Wo )[4/i(.g) + 2f 2 (g)} and J 2 = 
(J 2 M,)/i(ff) with h{g) = En=i9 2n /(nln) and f 2 (g) = 
£~=i Y.Z=i9 2(n+m) l[^-m\{n + m)]. In Fig. [TJ we plot 
the ratios J 2 / J\ and J z / J\ for various values of g and adi- 
abaticity parameter t/ujQ. We note that both (J 2 /Ji) x 
(u>o/t) and (J z /Ji) x (too/t) are functions of g only. 

The effective Hamiltonian in Eq. JB|) is different from 
that for spinless fermions in Ref. [2^. This is because 
the effective Hamiltonian for fermions contains an extra 
correlated hopping term J2j2jS8'^s^ n j c j+s' c 3+s (with 
Cj being the destruction operator for fermions) because 
the commutation relations for fermions are different from 
those of hcb given in Eq. ([3]). To see the difference 
clearly, let us consider the simplest case of one-dimension 
(ID). After carrying out the second-order perturbation 
theory for fermions (bosons), we get in ID the term 
c i j_ 1 CjCjCj + i (bj^bjbjbj+i) depicted by the process (a) 
in Fig. [5] and the term ctcj+ic'j^Cj (&t ftj+i &i ) de- 
picted by the process (b) in Fig. [5] For fermions, when 
these two terms are added, one gets cj_ 1 (l — 2n.j)cj+i 

whereas for bosons one gets only b^_ 1 bj + i. These ar- 
guments can easily be extended to d-dimensions. In a 
previous work 2 ^, while performing a similar second-order 
perturbation theory, the authors missed the process de- 
picted in Fig. Hfb). Here we would like to point out 
that, as mentioned in Ref. [24], the small parameter for 
our perturbation theory is t/(gu>o). 

III. MEAN FIELD ANALYSIS 

In this section, we shall study the phase transitions 
dictated by the effective Hamiltonian of Eq. © by em- 



3 




(b) Filled 



nian Hq is chosen as 

H = -J2^i-Si-B^2l, 



(10) 



where the molecular fields A* are obtained variationally 
by minimizing Fq. After some standard calculation, we 
obtain 

A* = A? = 2^Jx(5f + ,> + 2 £ J 2 (5^ 5 _ 5 ,> 
A* = 2B-J2j*( S i + s)o' ( n ) 



The eigenenergies of Eq. (fT0)l are 



where 



A = -£± A 



(12) 



FIG. 2: Depicted processes describe the following terms: (a) 
cj_ 1 Cjctcj+i when site j is empty; and (b) cjc^+ict ^ for a 
filled site j 



ploying the mean filed analysis (MFA) of Robaszkiewicz 
et a/.—. We first note that the hcb may be represented by 
spin one-half operators. More precisely, with the trans- 
formations S+ = S X + iS* = tf,S- =S X - iSv = b, and 
S z + 0.5 = b'b, the commutation relations of Eq. (|3|) are 
preserved. We can then write the Hamiltonian of Eq. Q 
in the following form: 



3,6 



— J2 2.^1 V a 3+6' a j+6 "+" a j+S' a j 

+o.5j z j2 s i s Us- B J2^ 2S 3 z + 1 )' ( ? ) 

i,<5 j 



with the constraint 



(8) 



where N is the number of sites in the lattice. Here, B = 
J z /2 + g 2 uo/2 is the effective magnetic field and n — 

iM.^i) i s the filling fraction (0 < n < 1). In the 
MFA, for a trial Hamiltonian Hq the following identity 
holds: 

F < F = ~]n Tr [exp(-0H o )] + (H - ff } , (9) 

where /3 = 1/ksT and (. . . ) is the thermal average with 
respect to the trial Hamiltonian Hq . The trial Hamilto- 



A' 1 .' 



(13) 



Here £ = a,b represents the two sub-lattices. The eigen- 
functions are given by 



,4- / f C \ i 1 ■ f U C 

W = cos ( y )lg > +sm ( ' 



1 



2 " 2 



> 



^ ( = -sm(y)|- > +cos(^)| - - 



>, (14) 



where sin6^ = 52"" and cosfl^ = j^-. At T = OK, ground 
state expectation value of S x and S z are given by 



and 



sin^r 



cosOr 



(15) 



(16) 



Now, to obtain the ground-state phase diagram we cal- 
culate the ground-state energy to be 

9 J[N 

= — sin^asin^f, ^_sin 2 # a ^2_sin 2 6>h + 



8j; 



8j; 



J B 

-^rcos6 a co$,6 b - —p- [cos6» Q + cos6» b + 2] , (17) 

where j[ = Z nn J 1 , j' 2 = Z mm J 2 , and j' z = Z nn J z , with 
Z nn and Z nnn being the number of nearest-neighbor (NN) 
and next-nearest-neighbor (NNN) hopping processes re- 
spectively. For the Hamiltonian of Eq. ((6|) , Z nn = 4 while 
Znnn = 12 because of the diagonal hoppings given by the 
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FIG. 3: Mean field phase diagram with Jl n = J„/2- 3J 2 . 

third term on the right-hand-side of Eq.®. Minimiza- 
tion of Eg, with respect to 9 a and Ob, gives the following 
two conditions: 

2£?sin6' Q = ^cos^sin^ + ^sin^aCOS^a + 0.5 J 2 sin6> a cos6>fc 
2Bsm6b — JjCos^sin^Q + J 2 sin0f,cos0f, + 0.5 J^sin^cos^a . 

(18) 

Now we note that, for a charge density wave (CDW) 
state we have (6 a ,&b) = (0,tt) or (6 a ,0b) — (tTjO); for a 
superfluid (SF) state a = Ob', and for a phase separated 
(PS) regime a ^ Ob and a ,b ^ or n. Then, from 
Eci. tfTg]) , we obtain the following expression for the phase 
boundary: 

A_J^A_3^ = 1 + (2n - l) 2 
2j; j; 2Ji Ji l-(2n-l) 2 ' 1 ' 

From Eq. fj 19[) . we see that we obtain the mean-field 
phase boundary of Ref. [25[ when J 2 = 0. We further 
note that NNN hopping does not change the qualitative 
feature of the phase diagram (see Fig. [3|); it only in- 
creases the critical value of J z / 1 J 1 at which the transition 
from SF state to PS or CDW state occurs. 



IV. DIAGONAL LONG RANGE ORDER AND 
STRUCTURE FACTOR 

Diagonal long range order (DLRO) is the typical prop- 
erty of a crystalline solid and information about the pe- 
riodicity in the solid is contained in the structure factor. 
For crystalline solids, the structure factor shows delta 
function peak at the reciprocal lattice points. In terms 
of the particle density operators the structure factor is 
given by 

5(q) =J2e^-^({n i n i ) (n,) (nj». (20) 

ij 

In this paper, we have calculated the structure factor 
using exact diagonalization technique for lattice clusters 
of size 4x4, Vl~8 x Vl8, and V20 x V20 . 
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FIG. 4: Normalized structure factor S""(7r,7r) = 
S(-K,n)/S max (TT,n) for three different lattice clusters at 
half-filling. Here, adiabaticity parameter t/uo = 0.1 and 
NNN hopping J 2 = 0. 
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FIG. 5: Depiction of S*(n,n) for three different lattice clus- 
ters at half-filling with t/u>o = 1.0 and J2 = . 



A. ,h = 

We will now present the structure factor results when 
J2 = 0, i.e., for the xxz- model. In Fig. [4] we 
have plotted the normalized structure factor S*(ir, n) — 
S{ir,n)/S max {Tr,ir) where S max (n, tt) corresponds to all 
particles in only one sub-lattice. The calculations were 
done at half-filling and for different lattice clusters with 
the adabiticity parameter t/u>o = 0.1. From Fig. 21 we 
see that, at half-filling, the system makes a transition to a 
CDW state at a critical b-p coupling strength g c w 2.15. 

On the other hand, at a larger t/u>o = 1-0, the transi- 
tion for a half-filled system occurs at a significantly lower 
value of g c » 0.9 (see Fig. [5]). This is because, for J2 = 0, 
the transition is governed only by the ratio J z / J\. Since 
{Jz/ Ji) x (wo/t) ' s a monotonically increasing function 
of g, for a larger value of t/u>o, it takes a lower value of 
g to attain the same value of J z /Jx- Another important 
point to note from Figs. [3] and \5\ is that S*(tt, tt) is al- 
most identical for different lattice clusters. The jump in 
the structure factor becomes sharper as we increase the 
system size. However this does not change the point of 
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FIG. 6: Comparison of S*(tt,tt) values for t/uia = O.l&l.O. FIG. 8: Plots of S*(tt,tt) for three different lattice clusters at 
Plots are for 5 particles in a 4 x 4 lattice and J 2 = . half- filling, t/oj = 1.0, and J 2 / . 




g FIG. 9: Comparison of S*(ir,ir) plots at two extreme values 

of t/ujQ. Figures are for 5 particles in a 4 x 4 lattice cluster 
FIG. 7: S*(n,n) for three different lattice clusters at half- an( j J 2 ^ . 
filling, extreme adiabaticity (t/wo = 0.1), and non-zero NNN 
hopping J 2 . 

B. J 2 



transition significantly. This means that the 4x4 lattice 
cluster is enough to have a reasonable estimate of the 
transition point. 

Next, we proceed to analyze the system away from 
half- filling. Without actually presenting the details of the 
calculations, we first note that, for N p < 4 in a 4 x 4 lat- 
tice, there is no evidence of a phase transition. Here, we 
present the results for total number of particles N p = 5 
in a 4 x 4 cluster. From Fig. [6l we see that the qualita- 
tive features of the transition are similar to those of the 
half-filled case. However, in detail, the two cases differ 
in the following sense. Firstly, the critical values of the 
b-p coupling are larger for N p — 5 with g c as 2.45 for 
t/tJo = 0.1, while g c w 1.70 for t/u>o — 1.0. Secondly, 
the 5-particle system never attains a fully CDW-state 
as seen from S*(tt,t:) being noticeably less than unity. 
Lastly, in Fig. [6j we see that S*(tt,tt) decreases slowly 
after attaining a peak value which is perhaps due to some 
special correlations which need to be examined. 



In this sub-section, we shall consider the effect of the 
additional NNN hopping J 2 . At half filling, for a small 
value of the adabiticity parameter t/u>o — 0.1, we find 
that the system undergoes a phase transition from a SF- 
state to a CDW-state at a critical boson-phonon coupling 
strength g c ss 2.16 (see Fig. [7]) which is very close to the 
case when J 2 = (see Fig. 0]). This is because, when 
t/ojQ is small, the ratio J2/J1 << 1 for all values of g 
(see Fig. [J). 

When we increase the value of t/u>o , as is evident from 
Fig. the value of the critical coupling g c decreases. 
The physical reason for this has already been discussed 
in section HV Al At half- filling, for t/uio — 1.0, we find 
g c ~ 1.53. Furthermore, for a given value of t/u>a, the 
system makes a transition to the CDW state at a lower 
value of g when NNN hopping J 2 = 0. This is in accor- 
dance with the mean field analysis, which shows that the 
presence of J2 delays the transition [see Eq. (fT^l) ]. This 
is because, for relevant values of g, when t/u>o = 1, the 
J2/J1 term is small but not negligible. Presence of J2 in- 
troduces disorder in the system. Hence, it takes a higher 
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J z /J\ ratio (i.e., a higher value of g) to make the system 
ordered. Similar to the J2 = case in a 4 x 4 lattice, we 
also find that for fillings up to 0.25 there is no transition 
to a CDW state while a CDW transition does occur for 
5 hcb. For N p = 5, as seen from Fig. \§\ g c ~ 2.45 for 
t/tJo = 0.1 while g c w 1.85 for t/uj = 1.0. 



V. OFF-DIAGONAL LONG RANGE ORDER 

The concept of off-diagonal long range order (ODLRO) 
was introduced by Penrose and Onsager- to understand 
the nature of the order in superfluids. Bose-Einstein con- 
densate is one example which shows ODLRO. Following 
Refs. [26| and [27|, we define the general one-particle 
density matrix as 



which can be re-written as 

/EiPM) £ip(*.2) 

«2,1) 0(2,2) 



V p(N,l) P(N,2) 



0(2, N) 



p(N,N) J 



(24) 



From Eq. ([24]) . it is easy to see that, for a translation- 
ally invariant system, (no) is indeed an eigenvalue of the 
one-particle generalized density matrix. We found in our 
calculations that, for all the relevant regions of the vari- 
ous parameter spaces, the value of (no) obtained accord- 
ing to Eq. ([2"2"| and the highest eigenvalue of the density 
matrix coincide quite accurately. 



B. Superfluid fraction 



P~(i,j) 



bib,. 



J_ V pCk'Ri-q-Rj) 

N ^ 

k,q 



KK 



(21) 



Here () denotes ensemble average and b[. is the creation 
operator for hcb in momentum space. It is easy to see 
that p becomes the diagonal one-particle density matrix 
when i = j . 



Condensate fraction 



It follows from Eq. ([2T[) that, for a translationally 
invariant system, 



(22) 



where no is the occupation number for the k = momen- 
tum state. Eq. ([2"Tj) gives the Bose-Einstein condensate 
fraction as 



n b 



E 



P~(i,j) 
NN„ ' 



(23) 



In general, to find rib, one constructs the generalized one- 
particle density matrix p and then diagonalizes it to find 
out the largest eigenvalue. This procedure alone does not 
tell us which momentum state corresponds to the largest 
eigenvalue. To find out whether the k = momentum 
state is macroscopically occupied or not, we proceed as 
follows. First let us see if (n ) is one of the eigenvalues. 
For this consider the following single-particle generalized 
density matrix, 



( 0(1,1) 0(1,2) 
0(2,1) 0(2,2) 



0(1,^0 \ 
0(2, N) 



\p(N,l) p(N,2) ••• p(N,N) J 



To characterize a superfluid, another important quan- 
tity of interest is the superfluid fraction n s . The order 
parameter for a superfluid is a complex number and it is 
taken to be (b) = yj (n ) /N e %e where the lattice constant 
has been taken to be unity. Spatial variation in the phase 
6 will increase the free energy (or simply the energy at 
T = 0) of the system. We consider an imposed phase 
variation that is a linear function of the phase angle, i.e., 
we take 8(x) = 60% where L is the linear dimension of 
the system along the x-direction. For simplicity we have 
chosen the variation in 9 to be only along the x-direction. 
With these considerations we can write the change in en- 
ergy to be 



E[0 O ] -E[9 = 0] = ^mN p n s 



-W(x) 
m 



(25) 



where E[6q] corresponds to an imposed phase variation 
when #0 7^ 0- Here, it is important to note that #0 should 
be small, because a larger 6q can induce other excitations 
which can destroy the collective motion of the superfluid 
component (see Ref. [28| and the references therein for 
details). We then get the superfluid fraction to be 



N 



N„t 



P l eff 



E[6 } - E[6 = 0] 



(26) 



where t e ff — h 2 /2m. For our Hamiltonian in Eq. 
we find t e ff = J\ + 8J2. Now, to introduce the phase 
variation, we impose twisted boundary conditions on the 
many-particle wave function. A twist in the boundary 
conditions is gauge-equivalent to modifying the hopping 
terms in the Hamiltonian of Eq. ^ . With this modifi- 
cation, the effective Hamiltonian becomes 



-g 2 u 



y^,,, ./ ; ^V"* " •7,;/,,.., 

3 3,S 

T i9£-(R ...1— R.j+s)i+ , 

J 2 2^ e 3+ b j+S' b 3+S 

0.5J z ^n,(l-n J+ 5), (27) 

3,S 
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t/O> =0.1 

J, = 



• 5 3.0 * T 



. • g 2.5 • 

• 5 2.4 

^-g< 2.3 x ' 

,.-tf = 2.2 
,-• ...» 

. g 2.1 




FIG. 10: Free energy plots at different fillings for J2 = 0, 
t/uio = 0.1, and various b-p couplings g . 



FIG. 11: Free energy versus particle number N p for parame- 
ters J2 = 0, t/uio = 1.0, and g at different values. 



where x is a unit vector in the ^-direction. 



VI. ENERGY CONSIDERATION 

In this section, we shall examine the possibility of 
phase separation for the system of hcb coupled with op- 
tical phonons. To this end, we have calculated the free 
energy for different number of particles in a 4 x 4 lattice. 
In Sec. |TV] for both J 2 = and J 2 > 0, we observed 
that the system is always a pure superfluid for N p < 4 
in a 4 x 4 lattice and that it is either a pure CDW or 
a pure superfluid at half-filling. After plotting the free 
energy at different fillings, if it is found that the curve is 
convex at a given filling, then the system at that filling is 
said to be stable; whereas, if the curve is concave at that 
filling, then the system would be unstable against phase 
separation. This procedure of calculating free energy at 
various fillings to figure out the stability of a system is 
equivalent to the well-known Maxwell construction. 



A. J 2 = 

In Sec. IIV Al we saw that, at half-filling, J 2 = 0, and 
t/ojQ = 0.1, the system makes a transition to a CDW 
state at g « 2.15. From Fig. [TQl we see that for g < 2.2, 
the system in between quarter-filling and half-filling is 
stable due to the convexity of the energy curve here. As 
g is increased, the system close to half-filling becomes un- 
stable first and then the lower-fillings becoming unstable 
progressively. For J 2 — 0, we see from Figs. [TUJ and 
[IT] that the qualitative nature of the free energy curves, 
with respect to phase separation, does not depend on 
the adiabiticity parameter. As expected from the expla- 
nation in Sec. HVi these figures show that the critical 
value g c (where the phase separation starts) decreases as 
t/u>o increases. We notice in Fig. [TU] (Fig. [TTj) that, 
for N p = 7, 6, &5, the energy curves become concave be- 
fore g = 2.3,2.4,&2.5 (g = 1.2, 1.4, &1. 7) respectively. 
Furthermore, for N p = 5, the phase separation seems to 



occur at approximately the same value of g at which the 
CDW transition occurs (see Fig. [6]). 



B. J 2 / 

For the extreme anti-adiabatic regime, the situation, 
when NNN hopping J 2 7^ 0, is not too different from the 
J 2 = case. This can be seen by comparing Figs. H"2"land 
[TOl drawn for t/u>o = 0.1. This is expected because, for 
the extreme anti-adiabatic regime, the ratio J2/J1 << 1. 
However, when t/uiQ — 1.0 (J2/J1 ratio is not negligible 
at values of g considered in Fig. I13p . the situation is quite 
different from the J 2 = case away from half-filling. In 
Fig. [UJ we saw that the structure factor revealed a CDW 
transition at g c « 1.85 for N p = 5. However, for J 2 7^ 0, 
the phase separation occurs at a higher value of g w 2.1 
as can be seen from Fig. 1131 The interesting implications 
of the CDW transition occurring before the PS transition 
will be discussed in the next section. 

One additional important feature for the J 2 ^ case, 
compared to the J 2 = case, is that the phase separation 
first occurs at the low-filling side. For 6 particles, the 
phase separation sets in at g f=a 3.0; while for 7 particles 
the corresponding g-value is expected to be even higher. 
We could not obtain the exact (/-value for PS instability 
for N p = 7 as our code had convergence problems when 
we tried to go beyond g — 3.0. 



VII. RESULTS AND DISCUSSIONS 

Here, we will analyze together, in one plot, the quanti- 
ties S*(tt, it), rib and n s that were presented in earlier sec- 
tions. For a half- filled system at J 2 = and t/wo — 0.1 
{t/ojQ — 1.0), we can see from Fig. [JJ] (Fig. \T5§ that 
the system undergoes a sharp transition to an insulating 
CDW state at g c w 2.15 (g c w 0.9). At g = g c , while 
there is a sharp rise in the structure factor S(ir, tt), there 
is also a concomitant sharp drop in both the condensation 
fraction nb and the superfluid fraction n s . Furthermore, 
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FIG. 12: Free energy at different fillings for t/ujo = 0.1 and 
different values of g. Here NNN hopping J 2 7^ 0. 
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FIG. 13: Free energy versus filling with t/uio = 1.0 and b-p 
coupling g at different values. Here too J2 > 0. 
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FIG. 15: Comparison of S*(n,n), rib, and n s for 8 particles 
when J2 = and t/bJo = 1.0. 
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while n s actually goes to zero, rib remains finite [as fol- 
lows from Eq. ([2"5]) ] at a value 1/N = 1/16 which is an 
artifact of the hniteness of the system. The lower critical 
value of g at higher values of the adiabaticity parame- 
ter t/ajQ has already been explained in Sec. HVl Thus 
at half-filling, in the absence of NNN hopping, a system 
of hard core bosons coupled with optical phonons under- 
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FIG. 14: Comparison of normalized structure factor S*(n,iv), 
condensate fraction rib, and superfluid fraction n s for 8 par- 
ticles when J2 — and t/uio = 0.1. 



FIG. 16: Comparative plots of S*(7t,7t), nj,, and n s for 8 
particles when t/uoo = 0.1, but J2 > 0. 



goes a transition from a superfluid state to an insulating 
CDW state. For a half-filled system, the presence of NNN 
hopping does not produce a qualitative difference in the 
plots, except for changing the critical value g c of tran- 
sition and that too only for adiabaticity values t/uj of 
the order of unity. This can be seen from Figs. [16] and 
1171 Here also rib and n s behavior complements that of 
S(tt, 7t); the values of rib and n a drop noticeably when 
S(ir, tt) increases sharply. These results for half-filling, 
with J 2 = and J 2 7^ 0, were already qualitatively pre- 
dicted in the mean-filed analysis of Sec. IIIII 

Away from half-filling, the system shows markedly 
different behavior compared to the half-filled situation. 
From Figs. US and [H for J 2 = 0, although S(n,Tr) dis- 
plays a CDW transition at a critical value g c , rib does not 
go to zero, again due to finite size effects, even at large 
values of g. For N p = 5 and t/uo = 0.1 (t/ujo — 1.0), 
we obtain the critical value g c = 2.45 (g c — 1.75). From 
Figs. [TO] and [TT] we see clearly that at these critical 
values of g, the free energy curves become concave for 
N p = 5. This suggests that the system is in a phase- 
separated state, i.e., it is an inhomogeneous mixture of 
CDW-state and superfluid-state. Thus away from half 
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FIG. 17: Comparative depiction of S*(n,ir), rib, and n s for 8 
particles when t/u>o = 1.0 and J2 > 0. 
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FIG. 18: Concomitant transitions depicted by S*(n,n), rit, 
and n s for 5 particles when J2 = and t/uo — 0.1. 

filling, when J2 = 0, our hcb-system undergoes a tran- 
sition from a superfluid-state to a phase-separated-state 
at a critical boson-phonon coupling strength. 

In the presence of NNN hopping and in the extreme 
anti-adiabatic limit also, the system's behavior for N p = 
5 is very similar to that of Ji — at the same adiabatic- 
ity as can be seen by comparing Fig. [20] with Fig. IT51 
and Fig. [TO] with Fig. [TO] However for t/uio not too 
small, when NNN hopping is present, the system shows 
a strikingly new behavior for a certain region of the g- 
parameter space. Let us consider the system at N p = 5, 
t/u} = 1.0, and J 2 7^ 0. Fig. |2"T1 shows that, above 
g w 1.85, the system enters a CDW state (as can be 
seen from the structure factor); however, it continues to 
have a superfluid character as reflected by the finite value 
of n s . Furthermore, Fig. [TO] reveals that the system is 
phase-separated only above g = 2.0. This simultaneous 
presence of DLRO and ODLRO, without any inhomo- 
geneity (for 1.85 < g < 2.1), implies that the system is a 
supersolid. Similarly, for 6 and 7 particles as well, we find 
that the system undergoes transition from a superfiuid- 
to a supersolid-state and then to a phase-separated-state. 
This is displayed in the phase diagram given in Fig. 1221 
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FIG. 19: Simultaneous phase transitions shown by S*(n,n), 
rib, and n s for 5 particles when J2 = and t/uio = 1-0. 
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FIG. 20: Comparative study of S*(n,n), nj,, and n 3 for 5 
particles when J2 7^ and t/uo — 0.1. 

Another point to be noted here is that, in Figs. [17] and 
[5T](i.e., for t/u = 1.0, J2 > 0, and small values of g), n s 
becomes smaller than rib- We feel that this is an artifact 
of the approximation used for the mass in Eq. (|26p . 

Finally, we shall present the interesting case of J\ = 
as a means of understanding the supersolid phase in the 
phase diagram of Fig. [55] The physical scenario, when 
Ji can be negligibly small compared to J2 (see Ref. [24| 
for a one-dimensional example), and the detailed results 
will be published later—. Here we only present the re- 
sults that are relevant to the conclusions made in the 
above discussions. It is quite natural that, when J\ = 0, 
all the particles will occupy a single sub-lattice for large 
values of nearest-neighbor repulsion. For a half-filled sys- 
tem, above a critical point, all the particles get localized, 
resulting in an insulating state. This can be seen from 
Fig. [23] One can see that (at J z /J2 ~ 7.2) the struc- 
ture factor dramatically jumps to its maximum value, 
while n s drops to zero and rib takes the limiting value of 
1/16 for reasons discussed earlier. This shows that above 
Jz/Jz — 7-2, the system is in a insulating state with one 
sub-lattice being completely full. However, away from 
half-filling, the system conducts perfectly while occupy- 
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FIG. 21: Comparison of S*(7r,7r), n;,, and n s for 5 particles 
when J2 7^ and t/wo = 1.0. 



FIG. 24: Comparison of 5'*(7r,7r), rii, and n 3 for 5 particles 
when t/uo = 1.0. Here too J\ — 0. 
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FIG. 22: Phase diagram depicting various phases for J2 > 
and t/uio — 1.0 



FIG. 25: Plot of free energy for different number of particles 
at various values of g when Ji = 0. 



ing a single sub-lattice because of the presence of holes 
in the sub-lattice. For instance, from Fig. [24] drawn 
for N p = 5, we see that the structure factor jumps to its 
maximum value at J z j Ji ~ 7.5, while n s drops to a finite 
value which remains constant above J z / J2 — 7.5. From 
Fig. [25] we see, based on the curvature of the free energy 
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FIG. 23: Comparison of 5*(7r, 7r), rib, and n s for 8 particles 
when t/uio = 1.0, but Ji = 0. 



curves, that the 5-particle system does not phase sepa- 
rate both above and below the transition. In fact, this 
single-phase-stability is true for any filling. This means 
that, at non-half filling and above a critical J z /J2, the 
system is homogeneous with simultaneous existence of 
both DLRO and ODLRO, i.e., the system exhibits su- 
persolidity! Thus, except for the pathological case of one 
particle, the system undergoes a first-order phase tran- 
sition from a superfluid to a supersolid state away from 
half-filling. 
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